0.1 Intro

En esta segunda clase aprovecharemos para aprender una serie de herramientas sumamente útiles al momento de trabajar con datos geográficos:

  • Geocodificar direcciones
  • Cargar capas base para dar contexto
  • Realizar mapas de calor o heatmaps
  • Reproyectar capas

0.2 Geocodificando direcciones

Muchas veces la información que tenemos hace referencia a un lugar en el terreno pero no tenemos un par de coordenadas para ubicar en un mapa. Por ejemplo si contamos con un listado que tiene direcciones, por ejemplo Tres Arroyos 614, que sabemos dónde está pero no tenemos un par de coordenadas que permiten ubicarla en un mapa.

Por suerte esto tiene una solución! Geocodificación es el proceso que permite transformar una dirección exacta en un par de coordenadas, permitiendo así tener un objeto geográfico con una geometría precisa.

Usar los servicios de Open Street Map, Google o Here son algunas opciones que tenemos pero hoy vamos a usar el servicio web de geocodificación del Gobierno de la Ciudad de Buenos Aires, el cual es libre y gratuito.

Hay un montón de funcionalidades que tiene, pero hoy vamos a usar 2 principales que son muy útiles: normalización y geocodificación. Es decir que en un simple paso vamos a “corregir” el campo de una dirección para pasarlo al nombre oficial y luego transformar esa dirección normalizada a un punto en el mapa!

La librería que nos facilita el uso de este servicio se llama RUMBA, la cual la instalamos con el clásico install.packages("RUMBA"), y los datos que vamos a usar son un listado de direcciones de farmacias en la Ciudad de Buenos Aires.

#install.packages("devtools")
#devtools::install_github("bitsandbricks/RUMBA")

library(RUMBA)
library(tidyverse)
library(sf)

farmacias <- read.csv("https://raw.githubusercontent.com/martoalalu/clase-geo-salud/master/data/farmacias.csv")

head(farmacias)
##   id           direccion
## 1  1      SANTANDER 5101
## 2  2   AV RIVADAVIA 8900
## 3  3 JURAMENTO y cabildo
## 4  4  AV CORRIENTES 5000
## 5  5   AV CABILDO y Pico
## 6  6   AV RIVADAVIA 6202

El dataframe tiene una columna que tiene la dirección. Qué bien! Pero no tiene lat-long… Qué mal! Podemos geolocalizarla igual! Qué bien!

Con la función mutate_USIG_geocode podemos agregar las columnas de lon y lat y así tener nuestro tan preciado par de coordenadas. Es bastante simple, en primer lugar le indicamos el dataframe (aca es farmacias) y luego le especificamos el campo que tiene la dirección entre comillas. Y listo, ya tenemos un dataframe con par de coordenadas. Magia!

farmacias_geo <- mutate_USIG_geocode(farmacias, "direccion")

head(farmacias_geo)
##   id           direccion                address_normalised       lon
## 1  1      SANTANDER 5101              SANTANDER 5101, CABA -58.47806
## 2  2   AV RIVADAVIA 8900          RIVADAVIA AV. 8900, CABA -58.48926
## 3  3 JURAMENTO y cabildo JURAMENTO AV. y CABILDO AV., CABA -58.45670
## 4  4  AV CORRIENTES 5000         CORRIENTES AV. 5000, CABA -58.43616
## 5  5   AV CABILDO y Pico          CABILDO AV. y PICO, CABA -58.47436
## 6  6   AV RIVADAVIA 6202          RIVADAVIA AV. 6202, CABA -58.45410
##         lat
## 1 -34.66620
## 2 -34.63574
## 3 -34.56204
## 4 -34.60083
## 5 -34.54076
## 6 -34.62574

Y ahora sí podemos pasarlo a objeto geográfico y luego al mapa!

farmacias_geo <- farmacias_geo %>% 
  filter(!is.na(lon), !is.na(lat)) %>% 
  st_as_sf(coords = c("lon", "lat"), crs = 4326)

farmacias_geo
## Simple feature collection with 48 features and 3 fields
## geometry type:  POINT
## dimension:      XY
## bbox:           xmin: -58.50898 ymin: -34.66739 xmax: -58.3647 ymax: -34.54076
## CRS:            EPSG:4326
## First 10 features:
##    id             direccion                  address_normalised
## 1   1        SANTANDER 5101                SANTANDER 5101, CABA
## 2   2     AV RIVADAVIA 8900            RIVADAVIA AV. 8900, CABA
## 3   3   JURAMENTO y cabildo   JURAMENTO AV. y CABILDO AV., CABA
## 4   4    AV CORRIENTES 5000           CORRIENTES AV. 5000, CABA
## 5   5     AV CABILDO y Pico            CABILDO AV. y PICO, CABA
## 6   6     AV RIVADAVIA 6202            RIVADAVIA AV. 6202, CABA
## 7   7           MORENO 2502                   MORENO 2502, CABA
## 8   8 LIBERTAD y Libertador LIBERTAD y DEL LIBERTADOR AV., CABA
## 9   9              PERU 402                      PERU 402, CABA
## 10 10         santa fe 3944             SANTA FE AV. 3944, CABA
##                       geometry
## 1   POINT (-58.47806 -34.6662)
## 2  POINT (-58.48926 -34.63574)
## 3   POINT (-58.4567 -34.56204)
## 4  POINT (-58.43616 -34.60083)
## 5  POINT (-58.47436 -34.54076)
## 6   POINT (-58.4541 -34.62574)
## 7  POINT (-58.40158 -34.61331)
## 8  POINT (-58.38284 -34.58823)
## 9  POINT (-58.37466 -34.61281)
## 10 POINT (-58.41905 -34.58314)

Genial. Ahora nuestro objeto es geográfico 100%. Por último, al mapa.

ggplot() +
  geom_sf(data=farmacias_geo)

Ok. Logramos mapear los puntos pero no tenemos una referencia clara, ni siquiere sabemos si todos están en la Ciudad de Buenos Aires. Como que el mapa está vacío…

0.3 Datos con contexto: capas base

Cuando manejamos datos geográficos el contexto es algo muy importante, tener una referencia sobre dónde están ubicados es clave ya que nos sitúa y al mismo tiempo comunica mejor. Una manera de hacer esto es agregandole una capa base, algo asi como el lienzo sobre el cual “pintamos” nuestros datos.

Para hacerlo vamos a usar las librerías ggmap y osmdata.

Básicamente lo primero que tenemos que hacer es indicarle a R cuál es el limite, el boundin gbox del mapa que queremos descargarnos. Usamos la función getbb y le decimos que es “Ciudad Autónoma de Buenos Aires” ya que es el nombre oficial. Si estuvieramos haciendo un mapa de por ejemplo Rosario le indicaríamos “Rosario, Santa Fe, Argentina”.

library(ggmap)
library(osmdata)

bbox <- getbb("Ciudad Autónoma de Buenos Aires,Argentina")

head(bbox)
##         min       max
## x -58.53145 -58.33512
## y -34.70564 -34.52655

El resultado es una matriz que tiene las X e Y mínimas y máximas que delimitan a la Ciudad de Buenos Aires. Ahora lo que tenemos que hacer es pedirle que nos tragia la capa base con get_stamenmap, la cual tiene una serie de parámetros a especificar:

  • En bboxle indicamos los límites específicos, en este caso lo tenemos guardado en un objeto asi que ponemos el nombre con el que lo guardamos (bbox).

  • En maptype le indicamos el tipo de mapa que queremos descargarnos. Hay capas base que tienen el terreno, otros solo las calles, algunos en blanco y negro. Ahora nos estamos descargando el toner-background, pero si quieren otro tipo pueden chequear la página de Stamen Maps.

  • En zoom le indicamos el nivel de detalle que queremos. A más zoom, más detalle, pero también más pesado el archivo que nos descarguemos…

#Descargamos el mapa
CABA <- get_stamenmap(bbox = bbox,
                      maptype = "terrain",zoom=12)

Luego para volver a mapear en vez de usar ggplot() vamos a usar ggmap() pero no se preocupen la sintaxis sigue siendo la misma!

Llamamos ala nueva función y entre parentesis le indicamos el nombre del objeto que tiene descargado el mapa que usaremos como capa base.

ggmap(CABA)

Perfecto! Lo bueno es que si no nos gusta este mapa podemos probar con otro muy fácilmente, por ejemplo con toner-background.

CABA <- get_stamenmap(bbox = bbox,
                      maptype = "toner-background",zoom=12)

ggmap(CABA)

Y ahora sólo resta agregar los datos como una capa más.

ggmap(CABA) +
  geom_sf(data=farmacias_geo, inherit.aes=FALSE)

Genial! Ahora sí tenemos una referencia más clara!

0.4 Más allá del choropleth

Hagamos otros tipos de mapas. Tener ubicaciones exactas, es decir objetos que son puntos, nos permiten hacer otros tipos de mapas como de calor (heatmap), de burbujas (bubblemap) o de densidad (densitymap), entre otros.

Una de las problemáticas del coronavirus es su elevada letalidad en personas mayores de 65 años, lo cual hace que los geriátricos sean focos especialmente críticos. Hagamos algunos mapas para ver cuál es su distribución en la Ciudad!

0.5 Reproyectando capas

Un dato clave al momento de manipular objetos geográficos es saber su Sistema de Coordenadas de Referencia (CRS) y Proyección. Muchas veces nos pasa que al manipular más de un objeto geo los mismos se encuentran en distintas CRS lo cual nos impedirá poder hacer las operaciones que querramos.

Veamos un ejemplo!

## Reading layer `perimetro' from data source `http://cdn.buenosaires.gob.ar/datosabiertos/datasets/perimetro/perimetro.geojson' using driver `GeoJSON'
## Simple feature collection with 1 feature and 3 fields
## geometry type:  POLYGON
## dimension:      XY
## bbox:           xmin: -58.53152 ymin: -34.70529 xmax: -58.33515 ymax: -34.52649
## CRS:            4326